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Abstract —We consider load scheduling on constrained 
continuous-time linear dynamical systems, such as automated 
irrigation and other distribution networks. The requested loads 
are rigid, i.e., the shapes cannot he changed. Hence, it is 
only possible to shift the order back-and-forth in time to 
arrive at a feasible schedule. We present a numerical algorithm 
based on using log-barrier functions to include the state 
constraints in the social cost function (i.e., an appropriate 
function of the scheduling delays). This algorithm requires a 
feasible Initialization. Further, in another algorithm, we treat 
the state constraints as soft constraints and heavily penalize the 
constraint violations. This algorithm can even be initiahzed at 
an infeasible point. The applicability of both these numerical 
algorithms is demonstrated on an automated irrigation network 
with two pools and six farms. 

I. Introduction 

Scheduling problems arise in a variety of contexts. A 
peculiar scheduling problem is studied in this paper. It 
involves a constrained dynamical system and the processing 
of request to apply load, with a hxed but shiftable prohle, 
on this system across time. The goal is to optimize a social 
measure of sensitivity to scheduling delay while satisfying 
hard constraints. This problem is motivated by an aspect of 
demand management in automated irrigation networks [1], 
[2], [3], and may arise in other areas. The main challenge 
associated with this problem relates to the rigidity of the load 
request, whereby the construction of a feasible schedule can 
only involve shifting requests back-and-forth in time. Relax¬ 
ation of the rigidity requirement can lead to a formulation 
as a (large) linear program [1]. 

The formulation of the rigid load scheduling problem 
here distinguishes itself in the following ways. By contrast 
with [2], [3], a dynamics relationship between the load 
and the constrained system states are modelled. In [2], [3], 
only static capacity constraints are considered. The load 
scheduling problem considered in [1] does include dynamics, 
however this is modelled in discrete time. By contrast 
a continuous-time setting is employed in this paper. The 
discrete time formulation in [1] leads to a mixed-integer 
program, which is difficult to solve [4]. The continuous¬ 
time formulation here, on the other hand, gives rise to two 
gradient based numerical algorithms. The first involve log- 
barrier functions and thus a feasible initial point. The other 
uses a soft encoding of the state constraints, with heavy 
penalty on constraint violation, which does not require a 
feasible initial point. Both algorithms lead to only locally 
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optimal solutions due to the non-convexity of the scheduling 
problem. 

The rest of the paper is organized as follows. We first for¬ 
mulate the problem in Section The numerical algorithms 
are presented in Sections |V] and Section In Section |V] 
the applicability of the developed algorithms is numerically 
studied on an automated irrigation network with two pools 
and six farms. Finally, Section VI concludes the paper. 


II. Problem Formulation 

Consider the continuous-time linear time-invariant dynam¬ 
ical system 


m 

x{t) = Ax{t) + Bu{t) + ^ EiWi{t), a:(0) = xq, (1) 

i=l 

where x{t) S K""’ is the state of the system, u{t) G 
is the control input (e.g., the water-level references in auto¬ 
mated irrigation networks), and Wi{t) G 1 < i < m, 

is a prohle-constrained (e.g. on-off) input signal representing 
the scheduled load on the system corresponding to the supply 
of resources to customer i (e.g., the flow of the supplied 
water to each farmer in an irrigation network). Throughout 
this paper, we assume that the control signal u{t) over the 
planning horizon [0,T] with T G M>o is a discrete-time 
signal passed through a zero-order hold, that is, Ui{t) = 
for all 1 < z < and all kA < f < (fc + 1)A with a given 
sampling time 0 < A < T. Although slightly conservative, 
this assumption allows us to work with finite-dimensional 
optimization problems instead of more complicated optimal 
control problems. For all integers 0 < k < K := [T/A]—1 
and 1 < z < Uu, we dehne = ei[step(f — 

kA) — step(t — {k+ l)A)],Vf G [0,T], where the mapping 
step : M —>■ {0,1} denotes the Heaviside step function, i.e., 
step(f) = 1 if t > 0 and step(() = 0 otherwise. Moreover, 
Ci G is the column-vector with all entries equal to zero 
except the z-th entry which is equal to one. Therefore, we 
get u(t) =uo + Yh=i ar,kii,k{t)yt £ [0,r], where 

uq G is the steady-state control input. 

The customers submit demands (wi(f))tgR, 1 < z < to. 
These demands are rigid (i.e., their shape cannot be changed). 
Hence, our decision variables are the delays that correspond 
to shifting the requested demand across the planning horizon, 
i.e., we select > 0 so that Wi(t) = Vi{t — r^) for each 
1 < z < TO. In doing so, the goal is to ensure that that 
the state of the network x{t) stays inside the feasible set 
X = {x G I Cx < d]. We can write this scheduling 



Cx{t)<dyt e [o,T], 

K n„ 


problem as 

m 

min y^hi{Ti), (2a) 

s.t. Tj < Ti < Ti, Vi e {1,... ,m} (2b) 

x{t) = Ax(t) + Buit) 

m 

+ '^EiVi{t - Ti),x{0) = xq, (2c) 

(2d) 
e [0,T], (2e) 

^ — Iz 1 

uo+M<'«(i)<uo+it,Vi G [0,T], (2f) 

where and are the bounds on the scheduling delay for 
demand i, u and u are the bounds on the control signal 
deviations u{t) — uq, and the continuously differentiable 
mapping : K —>■ R captures the sensitivity of customer 
i to the delay for scheduling its demand. Throughout the 
next section, we implicitly assume that T is long enough so 
that the optimization problem in (|^ becomes feasible with 
a constant nominal control input (i.e., if the demands are 
separated from each other “to some degree”, the state of the 
system stays feasible without any effort). This assumption is 
made to make sure that we can always find a feasible initial 
condition for the numerical algorithm, proposed in the next 
section, by simply separating the demands from each other. 
Towards the end of this paper, we present another approach 
for solving our scheduling problem that avoids requiring a 
feasible initial condition by treating the constraints on the 
state as soft constraints. 

III. Numerical Algorithm 

In this section, we present a numerical algorithm for 
solving Q by adding the state constraints in ( |2dl i to the 
cost function using log-barrier functions. Let us define 

xo{t) = e^p{At)xo + f exp(A(f —/3))i?uod/3, 

t 

xlkit) = [ exp(A(f-/3))B^*,fc(/3)d/3,Vi G 

•^0 VfcG iV}, 

Xi{t)= / exp{A{t-/3))E^v^{l3)dl3,\/i G 
Jo 

Since the underlying system in 0 is linear and time 
invariant, the solution of the ordinary differential equa¬ 
tion 0 can be written explicitly as x{t) = xo(t) + 


- a). Now, we 
rewrite the optimization problem in (|^ as 

771 

y^E(ri), (3a) 

1—1 K riu 

s.t. x{t) = Xoit)+ 


can 


min 


k—1 1=1 


This optimization problem is still difficult to solve as we 
have to check infinitely many constraints; see Let us 
use the notation Cj, 1 < j < p, to denote the rows of the 
matrix C G We add the state consttaints in ( |^ to 

the cost function using log-barrier functions. This transforms 
the optimization problem in Q to 


min^ _ J{{n)T=i,{M7=0Li) 


(4a) 


\nu \K 

WUi.k 

S.t. Ti < Ti < Ti, Vz G {1,... ,to}, (4b) 

Ui < oii^k < Ui,yi G {1,... ,n„}, 

V/cG iV}, (4c) 

where 


J((A)Ei,(Kfc)”=i)Li) 

771 P pT 

= elogf-C; 


xoit) 

+ X] X] + '^x7{t-Ti) 


K 


k—1i—1 


1=1 


+ dz]dt 


in which e G R>o is an appropriately selected parameter. 

Remark 1: With increasing e, the optimal solution is 
pushed further from the boundary of the feasible set. There¬ 
fore, to recover the optimal scheduling, we need to sequen¬ 
tially reduce e and employ the solution of each step as the 
initialization of the next step. This would result in a more 
numerically stable algorithm; see the log barrier methods 
in [5]. 

Lemma 1: ^ continuously 

differentiable function. Moreover, 

-Cz{Ax'l{t-Tl)+EtVt{t-Tt)) 
e „ \K ' ' 


E 


^ Jo -C^x{t] (Ti)™ {{ai^k)7=i)Zi) + 


Z = 1 


dapi 




= E 


CzX^At) 


'o -C^x(t] (Ti)™ 1 , ((ai,fe)”=i)f^i) + J 


^dt 


where 


xit; (Ti)™i,((a*,fc)”“i)fe=i) 

K Uu 771 

= Xo{t) + 

k—1i—1 1=1 

Proof: First note that 



(3b) 

d 

dn 


Cx{f) < d,yt g*To,t], 

(3c) 

/ elog(-C'^a;(t; (tOZi, ((ai.fe)”=i)Ei) + d^dt 

T^< Ti< Ti, Vz G {1,..., m}, 

(3d) 

Jo 

Ui < a^,k < Ui,yi G {1,... ,n„}, 
VfcG iV}. 

(3e) 


F -C^dx{f{Ti)^^^,{{ai^k)7yi)Zi)/dn 

J, -CMt; (a)“ 1, ((a..fc)£i)f=i) + 4 










where 


dn 


= - r.) 


= -x'i {t - n) 

= -{Ax^t - Ti) + ElVl{t - Tl)). 


Similarly, we have 
d 


da 


j/ Jo 


[ elog(-C^x(i; ((ai,fe)”=i)Li) + dz)dt 

Jo 


' -Czdxjt- {{ai,k)'l=i)k=i)/daj,t 

-Czx{t; (r,)™ 1, ((ai,fc)r=i)f=i) + dz 


where 


da 4 




-x{t\ {Ti)T=i, ((ai,fe)r=i)fc=i) = 


The rest of the proof follows from simple algebraic manip¬ 
ulations. ■ 

Now, we can use Algorithmic (overleaf) to recover a local 
solution of Q. We can select the step sizes and 
using backtracking line search algorithm [5, p.464] and ter¬ 
minate the algorithm whenever the improvements in the cost 
function becomes negligible. Unfortunately, this algorithm 
requires a feasible starting point (because the argument of the 
logarithmic functions cannot become negative). We remove 
this assumption in the next section by proposing a numerical 
procedure that treats the state constraints as soft constraints. 


IV. Soft Constraints on States 


In the previous section, we were required to find a feasible 
initialization to be able to run Algorithm [C Here, we take a 
different approach by solving the optimization problem 


no P 

min + y 




(5a) 

s.t. Tj < Ti < Ti, Vi G {1,..., m} (5b) 
x{t) = Ax(t) + Bu{t) 

m 

+ '^EiVi{t-Ti), (5c) 
x(0) = xo, (5d) 

K 

uit)=Y^^ai^kii^k{t),yt G [0,r],(5e) 

k X? — 1 

u<u{t)<u,ytG[0,T], (5f) 


where -d G K>o is an appropriately selected constant. In this 
problem, we may violate the constraints Cx{t) — d < 0, 
however, the term heavily penal¬ 

izes such violations. For small values of V, this term also 
penalizes the states being close to the boundary (of the 
feasible set), however, as we increase d, this term approaches 
zero inside the feasible set and infinity outside of the feasible 
set. 


Note that similar to the previous section, we can trans¬ 
form (|^ into 

P rT 


nv t' pi 

min Y'h^{Ti) + ^ 


(6a) 
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s.t. a;(f) = a;o(f)+ 


k—1i—1 


+ (6b) 

2—1 

Ui < ctt^k < G {1,.. ■,nu}, 

Vfc G a:}. (6c) 


Let us define 


xoit) 


-dz df. 


= Y] h,{T,) + V / exp ( V ( C; 

i=i z=ido V V 

K Uu m 

k—1i—1 i—1 

Hence, we may rewrite ^ as 

. ^ 0 ^ 4 K 'd'{{Ti)Zi,{{at,k)7yi)k=i), (7a) 

s.t. Tj < Ti < Ti, Vi G {1,... ,m}, (7b) 

Ui < cti,k < Ui,\/i G {1,.. .,nu}, 

VfcG TV}. (7c) 

Similarly, we can prove the following result regarding the 
augmented cost function. 

Lemma 2: J'{{Ti)'^i,{{ai^kYlZi)k=i) a continuously 
differentiable function. Moreover, 

^f{{r 4 )T=,A{a^,k)Zi)k=i) = 

- X / 'deypp{d{Czx{t] (Ti)™ 1, ((ai,fc)”=i)fci) - 4)) 

z=ido 

X Cz{Ax}(t - Tl) + Eivt{t - Ti))dt, 

^ J'{{r^)T=lA{o^^,k)Zl)k=l) 


da 


P pT 


= X / ^exp{d{C^x(t;{n)Zi,iia,^k)7=i)k=i) - dz)) 

z=ido 


X Czx'li{t)dt, 


where x{t; {{ai^k)7=i)k=i) defined as in 

Lemma [T] 

Proof: First, note that 

X exp{d{Czx{t; - dz))dt 

= f dexp{d{CzX{t; (n)™ i, {{a,,k)ifi)k=i) “ dz)) 

Jo 


X Cz 


dn 


x{t; {{a,^k)7fi)k=i) 


di. 










Algorithm 1 Projected gradient algorithm for scheduling rigid demands. 

Require; Feasible initialization (Ti[0])™;^ and 
1: for Z = 1 , 2 ,... do 
2: Update 


n\ 


_ 


Ml 1] A7 


(n)™! = (r.[/- l])£i 

_ 



V£e m}, 


and 




w =p: 


aj^e[l - 1] - 




i9a 






=1 = ((«i,fc[* - 1])” 


,Vj e 

V£e K}, 


where, for constants /3 < 7, = /? if a: < /3, P^[a;] = a; if /? < a; < 7, and P^[a:] = 7 if a; > 7. 

3: end for 


TABLE I 

Numerical parameters used in the simulation. 
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Fig. 1. The demand by the farms as requested (without the scheduling 
delays). 


Similarly, we have 

^ r expi^iCMt; {M7 =i)Li) - d.))dt 

pT 


= 7?exp(7?(C';,a:(f; (n)™!, 

Jo 

d 


X a 


da 






(Kfc)r=i)f=i 


)- 4 )) 

) df. 


The rest of the proof follows from simple algebraic manip¬ 
ulations. ■ 

Algorithm may be used with the gradients in LemmaJ^ 
to find a local solution of the optimization problem in (|^ 
Moreover if, after finding the optimal solution for a given J), 
the state constraints were violated at an intolerable level, we 
may sequentially increase and solve the problem until we 
get acceptable performance. 



Time ( min ) Time (min) 


Fig. 2. The shifted demands at the initialization of the algorithm. The 
dotted curve demonstrates the requests and the solid curve demonstrates 
their shifted counterpart. 


V. Numerical Example 

In this section, we illustrate the applicability of the algo¬ 
rithms on a water channel with two pools. The numerical 
example is borrowed from [1]. Each pool is modelled as 

yr{s) = -e q^{s) -g*+i s)-G s , 

s s s 

where Cin,i and Cout,i are discharge rates determined by the 

physical characteristics of the gates used to set the flow 

between neighbouring pools, and fd,i is the delay associated 

with the transport of water along the pool. Here, C,i{s) 

denotes the overall off-take flow load on pool i, that is, all the 

water supplied to the farms connected to this pool. Moreover, 

qi{s) is the flow of water from pool i — 1 to pool i and yi{s) 

denotes the water level in pool i. Eor the purpose of this 

example, we replace the delays with their first-order Fade 

approximatioif] Each pool is controlled, locally, by 




Ki{4>iS + 1) 

s{piS + 1) 


{ui{s) - y^{s)), 


'Note that the choice of a first-order Fade approximation is justifiable as 
the pool delays are all parts of closed-loops (with local controllers), with 
loop-gain cross-overs that are sufficiently small to make the overall closed- 
loop behaviour insensitive to the approximation error [6]. 
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Fig. 3. The shifted demands for the local solution recovered by 
Algorithm with e = 0.1. The dotted curve demonstrates the shifted 
demands at the initialization and the solid curve demonstrates the shifted 
demands at the locally optimal solution. 



Fig. 4. The reference signal for the local solution recovered by Algorithm^ 
with e = 0.1. The red lines show the boundary of the feasible region. 


where (j)i, and pi are appropriately selected control 
parameters. Furthermore, Ui{s) denotes the water-level ref¬ 
erence signal of pool i. Table [I] shows the parameters used 
in this example. The state constraints are as follows 9.4 < 
yiit) < 9.7 and 9.5 < 2 / 2 (f) < 9.7. Finally, throughout this 
example, we fix uq = [9.50,9.55]^. 

Figure [T] illustrates the requested demands of the farms. 
Here, and respectively, denote de¬ 

mands for pool 1 and 2. Let us select linear penalty functions 
hi{Ti) = Ti for all i. Moreover, assume that the reference 
signal should belong to a bounded region captured by 


Uq - 


0.05 

0.05 


< u{t) < Ug -f 


0.05 

0.05 


Note that without these control input constraints, one can 
schedule all the loads without any delay but with large 
control input deviations. In this example, we select = 0, 
Vt, which means that we can only shift demands forward. 
Let us also select Ti = 300 min for all i. 

First, we use Algorithm [T] to extract a reasonable schedule 
by shifting these demands. This algorithm requires a feasible 
starting point, which can be constructed by shifting the 


Fig. 5. The outputs for the local solution recovered by Algorithm^ with 
t = 0.1. The red lines show the boundary of the feasible region. 
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Fig. 6. The output of the system when all the decision variables (the 
control inputs and scheduling delays) are set equal to zero. 
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Fig. 7. The shifted demands for the local solution recovered by the 
proposed algorithm in Section |IV| with = 100. The dotted curve 
demonstrates the shifted demands at the initialization and the solid curve 
demonstrates the shifted demands at the suboptimal solution. 
































































































Fig. 8. The outputs for th e loc ally optimal solution recovered by the 
proposed algorithm in Section EYl with I? = 100. The red lines show the 
boundary of the feasible region. 


Fig. 10. The outputs for th e loc ally optimal solution recovered by the 
proposed algorithm in Section m with 'd = 10. The red lines show the 
boundary of the feasible region. 



Fig. 9. The reference signal for t he l ocally optimal solution recovered by 
the proposed algorithm in Section with 'd = 100. The red lines show 
the boundary of the feasible region. 


demands (to be somewhat distant from each other). Figure]^ 
illustrates the shifted demands at this initialization (solid 
curve) as well as the original requests (dotted curves) for 
comparison. Let us fix e = 0.1. Figure shows the shifted 
demands for the local solution recovered by the proposed 
algorithm. As we expect, the delays are significantly smaller 
in comparison to the initialization. Figure portrays the 
reference signals for the solution and Figure [^illustrates the 
outputs. Evidently, the output stays in the desired region. 
Although powerful. Algorithm requires a feasible initial 
condition that may not be easy to find. Therefore, in the rest 


of this section, we study the method presented in Section IV 


Let us select d = 100. Figure [^ illustrates the shifted 
demands for the locally optimal solution recovered by the 
proposed algorithm. At the starting point (fed to the algo¬ 
rithm), all the decision variables are selected to be equal 
to zero for which the state of the system does not stay in 
the feasible region; see Figure [^ Figures and [^ show the 
output and the control for the local solution recovered by 
the proposed algorithm in Section [W] Interestingly, all the 
constraints are satisfied, which is because of the large value 
of 1 ?. If we reduce to be equal to 10, the outputs violate the 
constraints on the state as shown in Figure Evidently, by 


comparing Eiguresj^j^ and[T^ we can see that by increasing 
■d, the constraint violations are becoming more infrequent 
(until they do not occur at all). 

VI. Conclusions 

In this paper, we presented numerical algorithms for 
scheduling demands on continuous-time linear time-invariant 
systems. The rigidity of the demands dictated that we can 
only shift them back-and-forth in time (and cannot change 
their shapes). The first algorithm used log-barrier functions 
to include the state constraints in the cost function. The 
second algorithm considered the state constraints as soft 
constraints and added a penalty function for the constraint 
violations to the cost function. Euture research can focus on 
constructing a market mechanism for achieving the optimal 
schedule based on the customers preferences. We can also 
compute an optimality gap via finding a lower-bound for the 
solution of the dual of the problem corresponding to 
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